Statistical hadronization with exclusive channels in e + e annihilation 
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We perform a systematic analysis of exclusive hadronic channels in e + e~ collisions at centre-of- 
mass energies between 2.1 and 2.6 GeV within the statistical hadronization model. Because of the 
low multiplicities involved, calculations have been carried out in the full microcanonical ensemble, 
including conservation of energy-momentum, angular momentum, parity, isospin, and all relevant 
charges. We show that the data is in an overall good agreement with the model for an energy 
density of about 0.5 GeV/fm 3 and an extra strangeness suppression parameter 7s ~ 0.7, essentially 
the same values found with fits to inclusive multiplicities at higher energy. 
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The statistical approach to multi-hadron production e + e _ annihilations has a long story. Early works date back 
to the '70s [1, 2], with different versions of the model and different observables examined in the relevant analyses, 
such as inclusive yields, multiplicity distributions etc. All of these calculations involved simplifying assumptions and 
drastic approximations, mostly because of the lack of computing power, so that in practice it was very difficult to 
confirm or rule out the statistical model, also in view of it being conceived as a full alternative to a dynamical model. 
Q_T Nowadays, with QCD being the accepted theory of strong interactions, the Statistical Hadronization Model (SHM) 
has resurged as a model of hadronization, it has a framework based on quantum statistical mechanics [3] , and it has 
been extensively and succesfully applied to the analysis of inclusive multiplicities in elementary [4] and relativistic 
heavy ion collisions [5]. Also, the SHM was shown to succesfully reproduce transverse momentum spectra in hadronic 
collisions [6] with specific predictions concerning the approximately exponential shape of the 1ow-pt spectrum and 
the so-called ttit scaling phenomenon. This model has by now become a standard tool in heavy ion physics while the 

. reasons of its success in elementary collisions are still subject of debate [7] . 

It is a common belief that statistical equilibrium cannot be attained in elementary collisions via post-hadronization 
collisions because of the low multiplicities and the rapid expansion. Therefore, one is led to conclude that statistical 
equilibrium is an intrinsic feature of the hadronization process itself, as envisaged by Hagedorn many years ago 

\ ("hadrons are born at equilibrium" [8]). In the latter case, two possibilities arise: 

Q\ 

• - the apparent statistical equilibrium is just mimicked by a special property of the dynamics governing the 
hadronization tending to evenly populates all final states, but which has essentially nothing to do with a proper 
statistical system, which can be realized only within a finite volume. 

• - the apparent statistical equilibrium is established within a finite volume and therefore is a "genuine" one; if 
. the volume was large enough, a proper temperature could be introduced. 

5_i 

The former picture can be defined as phase space dominance to discriminate it from proper statistical equilibrium. 
While phase space dominance is, though, a highly non-trivial hypothesis, its predictions do quantitatively differ from 
the proper SHM. As pointed out in refs. [3, 9], the dynamical matrix element of the decay of a massive cluster into 
N particles contain, in the SHM case, peculiar quantum statistics terms (Bosc-Einstein and Fermi-Dirac correlations) 
owing to the finite cluster volume, which are generally absent in the phase space dominance picture. The very fact 
that Bose-Einstein and Fermi-Dirac correlations have been observed in elementary collisions demonstrates the finite 
extension of the hadron emitting source and therefore favours a model like the SHM where finite volume is a built-in 
feature. 

Anyhow, it would be desirable to quantitatively test the the genuine statistical hadronization model on observables 
which arc more sensitive to the form of the matrix clement than average multiplicities. For this purpose, in this work 
we compare the production rates of exclusive channels in e + e~ collisions at low energy with the predictions of the 
SHM. 

Exclusive channels in e + e~ collisions have been measured at low centre-of-mass energy (< 4 GeV). At such a low 
energy, QCD is in the full non-perturbative regime and one can assume that, unlike at higher energy where clusters are 
two or more (jettiness), only one hadronizing massive cluster at rest in centre-of-mass frame of the collision is forme. 
The price to be paid is that, in calculating the model predictions, none of the relevant conservation laws, including 
energy-momentum, intrinsic angular-momentum and parity, as well as internal symmetries, can be neglected (see e.g. 
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ref. [10] where pp annihilation at rest was studied in the SHM). In the SHM framework, this means that one has to 
calculate averages in the most general microcanonical ensemble of the hadron-resonance gas. 

To carry out this calculation, in this work we take advantage of the formalism developed in two previous papers 
of ours [11, 12] where the microcanonical partition function of an ideal multi-species relativistic gas was calculated 
enforcing the conservation of the maximal set of observables pertaining to space-time symmetries (energy-momentum, 
spin, helicity, parity). We extend the formulae obtained therein to the hadron-resonance gas including internal 
quantities conserved by strong interaction (isospin, C-parity and abelian charges). We then take into account resonance 
decays and compare the results of our calculations with the data collected in e + e~ collisions at low energy. 

The paper is organized as follows: in Sect. II we will expound a formulation of the SHM in the full microcanonical 
ensemble which is suitable for the problem of exclusive channels and in Sect. Ill we will obtain an expression for the 
their rates; Sect. IV will be focussed on a method to compute them numerically. Finally, Sect. V will describe the 
analysis of data in e + e _ collisions at low energy will and Sect. VI will be devoted to a discussion of the results and 
to conclusions. 



II. THE STATISTICAL HADRONIZATION MODEL IN THE MICROCANONICAL ENSEMBLE 

In the modern formulation of the SHM [3] , the strong interaction process in a collision between particles leads to 
the formation of a set of extended massive objects called clusters or fireballs). Each cluster decays into hadrons in 
a purely statistical fashion, that is any multi-hadronic state within the cluster compatible with its quantum numbers 
is equally likely. The number of clusters produced, as well as their kinematical and internal quantum properties, 
are determined by the prior dynamical process and are not predictable within the SHM itself. Particularly, in high 
energy collisions (y^s > 10 GeV), the production of clusters following the perturbative parton shower stage leads to 
a multiple cluster production. Conversely, for energies sufficiently below the perturbative regime (y/s < 4 GeV), one 
may expect that, to a very good approximation, a single cluster is formed (see fig. 1). Under this circumstance the 
whole centre-of-mass energy is spent to produce particles and no jet is observed. The cluster mass will then coincide 
with yfs and its relevant quantum numbers will be the same as those pertaining to the initial state of the collision. 



FIG. 1. Dominant contribution to hadroproduction in a low energy e + e collision with formation of a single hadronizing cluster. 

Once the quantum numbers of a cluster are known, the assumption of equal probability allows to perform cal- 
culations within the framework of statistical mechanics, in the relativistic microcanonical ensemble. The use of the 
microcanonical ensemble is necessary as at such a low energy the effect of exact conservation laws is very important [13] . 
It should be emphasized that the basic assumption of the model states that multi-particle localized states compatible 
with cluster's conserved quantities are equally likely, but these states do not coincide with observable free particle 
asymptotic states. Such a distinction is, for practical purposes, not an issue when the volume is sufficiently large, 
but it is relevant in principle and may result in quantitative differences when the volume is small, i.e. comparable 
with the third power of the typical lenght scale of the hadron world, the pion Compton wavelength 0(1) fm 3 . This 
is discussed in detail in refs. [3, 9]. 

If the cluster can be described as a mixture of states, the basic postulate implies that the corresponding density 
operator is a sum over all localized states projected onto the initial cluster's quantum numbers through a projection 
operator P^: 





(1) 
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where \hv) are multi-hadronic localized states and Pj is the projector onto the cluster's initial conserved quantities: 
energy-momentum, intrinsic angular momentum and its third component, parity and the generators of inner symme- 
tries of strong interactions . If, on the other hand, the cluster is prepared in a pure quantum state |^>)(what is the 
case for a single produced cluster in e + e _ collisions) then, according to the basic assumption, this is ought to be an 
even superposition of all multi-particle localized states with the initial conserved quantities, that is: 

|V>) = y^Ch v Pj\hy) with \c hv \ 2 = const. (2) 

hv 

It can be readily shown [3] that if the coefficients Ch v have random phases, an effective mixture description with the 
operator in (1) is recovered. Hence, a new hypothesis is introduced: if the cluster is a pure state, the superposition 
of multi-hadronic localized states must have random phases. 

The operator Pj can be formally defined as the pseudo- projector (it is not idempotent, see below) onto an irreducible 
vector of the full symmetry group and worked out in a group theory framework [9, 12, 14]. It can be factorized into 
a "kinematic" pseudo-projector, associated to general space-time symmetries, and an actual idempotent projector 
for inner symmetries, associated to compact groups. For the general space-time symmetry the relevant group is 
the extended orthochronous Poincare group 10(1, 3) 1 " and an irreducible state with positive mass is identified by a 
four- momentum P, a spin J and its third component A and a discrete parity quantum number w = ±1. Therefore: 

Pi — Pp, J,X,7rPinner (3) 

If the pseudo-projector Pp,j.\ is worked out in the cluster's rest frame where P — (M, 0), it further factorizes [12, 14], 
i.e.: 

Pp,y,A, w = S 4 (P - P)Pj,xP* = S 4 (P - P)Pj,x-^- (4) 

where P is the four-momentum operator, Pj.\ is a projector onto SU(2) irreducible states \ J, A) and fi is the space 
reflection operator. Thus, the pseudo-projector (3) becomes: 

P, = 5\M - P°)5 3 (P)P J , A P 7r P mner (5) 

Note that Pj.x, P n and Pinner commute with each other. 

As clusters are colour singlets by definition, the projector Pinner involves flavour and baryon number conservation. 
In principle, the largest symmetry group one should consider is SU(3) flavour, plus three other U(l) groups for baryon 
number, charm and beauty conservation. However, SU(3) symmetry is badly broken by the mass difference between 
strange and up, down quarks, so it is customary to take a reduced SU(2)<g>U(l) where SU(2) is associated with isospin 
and U(l) with strangeness. The isospin SU(2) symmetry is explicitly broken as well, but the breaking term is small 
and can generally be neglected. However, most calculations in the past have replaced isospin SU(2) with another 
U(l) group for electric charge, so that the symmetry scheme, from an original SU(2) isospin (g)U(l) s t r angeness®U(l)baryon 
reduces to U(l) c hargc<8>U(l) st rangeness^U(l)b ar y n- 

Altogether, Pinner can be written as 

Pinner = PiJ^Pq^X (6) 

where / and I 3 are isospin and its third component, Q = (Qi, . . . , Qm) is a vector of M integer abelian charges 
(baryon number, strangeness, etc.) and P x is the projector onto C-parity, which makes sense only if the system is 
completely neutral, i.e. I — and Q = 0; in this case, P x commutes with all other projectors. 

From the density operator (1) the probability of observing an asymptotic multiparticle state |/) ensues: 

= mpvm (7) 

Pf tr(P 4 P v P 4 ) (<> 

which is a well-defined one with regard to positivity and conservation laws because pf = if the state |/) has not 
the same quantum numbers as the initial state. The normalizing factor, i.e. the trace of the operator PiPyP^ can be 
worked out as: 

tv(PiP v Pi) - tr(p2p v ) - «5 4 (0)tr(P l P y ) . (8) 



1 Operators in the Hilbert space will be denoted with a hat. Exceptions to this rule are projectors, which will be written in serif font, i.e. 
P. 
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where we have used the particular form of Pi in eq. (5) and we have taken into account that all operators except 
6 4 (P — P) are idempotent. The reason for the presence of a divergent positive constant S 4 (0) is the non-compactness 
of the Poincare group, which makes in fact impossible to have a properly normalized projector. The last trace in (8) 
can be written as 

tr(P i P v )=J2(hv\Pi\h v ) = n (9) 

hv 

which is, by definition the microcanonical partition function [11], i.e. the sum over all localized states projected onto 
the conserved quantities defined by the selected initial state. If only energy and momentum conservation is enforced, 
il takes on a more familiar form: 

n = Y,(hv\5 4 (P-P)\h v ). (10) 

hv 

In principle, the asymptotic multi-particle states |/) in eq. (7) only include strongly stable hadrons, while the inter- 
action between them is understood in the same equation through the projector Pj which contains the full hamiltonian 
P° (see eq. (5). An outstanding theorem by Dashen, Ma and Bernstein [15] asserts that, in the thermodynamic 
limit, the partition function - in any ensemble - of an interacting system is the sum of the partition function of the 
system without interaction and a term involving the scattering matrix between the otherwise free particles. The 
well-known consequence of this theorem is the so-called hadron-resonance gas model; if only the resonant part of the 
scattering matrix is retained (the background interaction is neglected), the interaction term of the partition function 
reduces to that of a gas of resonances treated as free particles with distributed mass. Strictly speaking, there is an 
additional contribution from resonance interference, which might be sizeable in case of wide, overlapping resonances 
with the same decay channel, but this depends on mostly unknown complex parameters and is thus assumed to vanish 
altogether or it is simply disregarded. 

In the spirit of the DMB theorem and the hadron-resonance gas model, we will therefore calculate the probabilities 
(7) including resonances as free particles with distributed mass in the multi-particle free states |/) and let them decay 
afterwards. It must be stressed that this is an assumption going beyond the scope of validity of the DMB theorem, 
which affirms the equality of two traces, and not of single trace terms. In other words, the hadron-resonance gas 
decomposition, strictly speaking, applies only to fully inclusive quantities and not to partly inclusive like multiplicities 
of single species or exclusive final states. Furthermore, the DMB theorem requires the thermodynamic limit V — > oo. 
Up to now, these problems have been ignored and the hadron-resonance gas model has been used to calculate hadron 
abundances and spectra when applying the SHM to the data. As has been mentioned, we will continue to use the 
hadron-resonance gas model in its simplest form also for exclusive channel rates and for small clusters. This is likely 
to be a good approximation but it should be kept in mind that deviations may well be implied. 



III. THE MICROCANONICAL CHANNEL WEIGHT 



Denoting a multi-particle final state |/) as \{Nj}, {p}) where {Nj} — (N\, . . . ,Nk) is the set of multiplicities Nj 
for each particle species j, i.e. the channel, and {p} stands for the set of kinematic variables (momenta and spin 
components or helicities) of the particles, we can calculate the probability of the channel {Nj} as: 

M 

= ^ip 154 ^ D^'}' MlP^P.P-erPyPiliJVil, {P}) (11) 

where use has been made of the decomposition in eq. (5) and the fact that the final states are eigenstates of total 
energy-momentum. One of the difficulties of working out this expression in a relativistic quantum field framework 
is that the localized states, in general, are not states with a definite number of asymptotic particles, unlike in non- 
relativistic quantum mechanics (see discussion in rcf. [11]). However, since we sum over all kinematic states and the 
projectors Pj.a, Pit do not change the number of particles, i.e. the set {Nj}, and we can use the ciclicity for these 
two projectors and rewrite the last expression as: 

= h^ i{Njl WI P — PvPiP.Pj,x\{Nj},{p}) (12) 
M 
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Now Pinner commutes with Py because localization does not affect internal symmetries, as well as with Pi, so we get: 

(13) 



P{N j} = ^({"ih MlPvPiPinner P.Pj.aK^}. W> = ^ £<W}, WIP^K^}, W) 



where we have used, in the last equality, the factorization (5) and the idempotency of projectors Pj nner , P n , Pj,\. 
Finally, the operators Py and Pi can be swapped in position using again the factorization (5), the commutation 
between Pinner and all other projectors, the ciclicity of P-k,Pj.\, and the fact that |{^Vj},{p}) are eigenstates of 
total energy- momentum. Hence, the relative probability of a channel {Nj} is proportional its microcanonical weight, 
defined as: 



P{Nj} cx n {Nj} = MIP.PvKiv,}, M) 

{p} 



(14) 



The microcanonical channel weight ^{a^} has been calculated explicitely in ref. [12] for an ideal relativistic gas 
of particles with spin with the full Poincare projector (4) in a quantum field theoretical framework. The obtained 
expression is essentially the same as the one would get in a multiparticle approach, i.e. working in the multiparticle 
tensor space with symmetrization for bosons and antisymmetrization for fermions; the only relevant quantum field 
effect being is an overall immaterial factor (0|Py|0). Let N be the total number of particles in the channel, i.e. 



^2f = i Nj — N; Sj and r\j respectively the spin and the intrinsic parity of the j-th particle species, p n the four- 
momentum of the n-th particle. Then, for a spherical cluster, the microcanonical channel weight reads [12]: 2 



1 hnAPj)' 



x n n F v\p Pj (nj) R^wpn^+mif n n 4*w»,) + ^^w^) 

where 



(15) 




n, 



K 

n< 



JV, 



(16) 



and p — (pi, . . . , pk) is a set of permutations, pj belonging to the permutation group x(Pj) is the parity of the 
j-th permutation and bj = 0, 1 if the species j is a boson or a fermion respectively; the symbol h n . (pj) in (15) stands 

for the number of cyclic permutation with rij elements in pj so that X^=i n j^ l n j (pj) 
Fourier integrals over a spherically symmetric volume, which for a sharp sphere read: 



N 4 



In eq. (15), Fy's are 



F^>( Pp{n) - R^WPn) = 



1 



(2tt) 3 



d 3 x e ix '( p p( n >~ R 3 1 (^)P") 



(17) 



= R 2 jidPpW-Rg'WP^I^) 

2 ^ 2 |Pp(n) - R3 'WPnl 

R being the radius, j\ the spherical Bcssel function of the first kind and P-3(ip) is a rotation of an angle tp along the 
z axis. The factor (0|Py|0), as has been mentioned, is immaterial as it cancels out in the ratios between different 
channels. 

In the cq. (15) the dependence on the cluster polarization state A has disappeared because of spherical simmetry 
[12]. 



2 We take this opportunity to notice that in the formula (82) in ref. 4 there was a factor 1/2 in excess. 

3 The set of integers hi, ... , /ijy = {hn}, is usually defined as a partition of the integer N in the multiplicity representation. 
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If, in eq. (15), we sum up over all angular momenta J and neglect all permutations except the identity (corresponding 
to the Boltzmann statistics), we obtain the more familiar expression: 



V 



N 



{N 3 } 



(2tt)37V 



(2S 3 + l)f 

3=1 



d 3 pi 



d 3 p w 5 4 



' N 



(18) 



n=l 



which can be used to show that the dynamical matrix element in the cluster's decay, according to the SHM, is 
proportional to P ■ pi/rho for each particle, p is the proper energy density of the cluster [3]. 



A. Internal symmetries 

The eq. (15) only contains the Poincare group projector (4); we now have to include the internal symmetry pro- 
jector Pinner- First of all, we will disregard altogether the projector Pq on the abelian charges baryon number and 
strangeness, as this gives rise to non-trivial coefficients (i.e. or 1) and can be easily implemented algorithmically 
just by enforcing y\ QjNj = Q, where q 3 is the vector of abelian charges for the species j. 

As has been mentioned, we will work in the multiparticle tensor space instead of Fock space and, for this purpose, 
it is convenient to introduce the concept of particle type. Particles of the same type are to be taken as identical, 
yet in a different charge state. We will take as identical (hence of the same type) light-flavoured non-strange mesons 
belonging to the same isospin multiplet or if they are a particle- antiparticle pair. For instance, 7r + ,7r~, and ir° belong 
to the same type, that one can define as the pion. Similarly, p and p or K + and K~ belong to the same type, while 
p and n, or K + and K - do not because, albeit forming an isodoublet, are not light-flavoured mesons. 

Let us denote with Lj the number of particles of the type j in a given channel {Nj} (a channel {Nj} completely 
defines the corresponding set {Lj}, while the converse is not true). If pj e Si is a permutation of the integers 
l,...,Lj, x(Pj) its parity and bj = or b 3 ■ = 1 if particles of type j are bosons or fermions respectively, the general 
final state |/) in the multiparticle tensor space can be written factorizing groups of particles of the same type as: 

K 

n {k Pjih) }, {v Pjih) }, k,( W }> . (19) 

3=1 

K being now the total number of types; p the set of permutations p\, . . . , pk, the kinematical variables of the 

particle lj (momentum and polarization), r^. its parity and q/^ its quantum numbers; I 3 3 is the isospin third component 
of the lj-th particle of the type j. If we define: 

\f p ) ee J] {qp i( ,,)}> (20) 

3 = 1 

then |/) reads, according to (19): 

I/p) (21) 

and the state weight ujj ee (,/|PiPy|/) for a channel {Nj} in the multiparticle tensor space can then be written in 
terms of the corresponding set {Lj} as: 

(A|PiPv|/p> • (22) 

where i is the identical permutation. 

We shall start studying the action of projectors on states with definite C-parity Xi I and ^3 bearing in mind that 
relevant projectors can be moved to the right of Py to act on | /) . As already stated, the projector Pc is meaningful 
only if the cluster is completely neutral and reads: 

Pc|/> = ^%) (23) 



i/) = E 



K 

n 

3=1 



XiPj 



fL-\ 



if) - E 



g x{p 3 ) b 



3 = 1 



EE 



K 

n 

3=1 



X{Pj°. 



'3 U 3, 

Lj\ 



Cf CT |p,Pvi./p) = E 



K 



\{x(P3) b 



3 = 1 
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where C is the charge-conjugation operator transforming the state |/) into: 



ci/> = E 



v L 3 l 

3 = 1 v J 



C\f P )= X cJ2 



n 



i/p) 



(24) 



In the above equation: 



K 



I/p> = LI ^AkMi^A-irndv^d-^h)}) 

3 = 1 



(25) 



and xc is the product of intrinsic C-parities of the completely neutral mesons and the charge conjugation phase factors 
of non-strange charged light-flavoured mesons, defined in Appendix A. The symbol rj l . stands for the parity of the 
charge-conjugated lj-ih particle. In order to ensure the commutation between Pc and the space reflection operator 
n, which has already been used in this Section, we have to set rj = rj for all particles. This is obvious for bosons but 
not for fermions, whose parity are arbitrary provided that \r]\ = 1 and rjrj = — 1. We then set TfB =rj^ = i for baryons 
in order to meet all requirements. 

Let us now move to the isospin projector which can be written formally as: 



P Ih = \I,h)(I,h\ • 



(26) 



Preliminarly, it is useful to write the state \f p ) as the tensor product of a ket \s) including kinematical variables k 
and parities 77, and an internal part: 



3=1 3=1 



(27) 



so that: 



i/> - E 



n 



xipj 



3=1 v j 



K 



\s P )j\\i3Ai p r s) }M Pj{h )}) 
3=1 



(28) 



Now, by using eqs. (23) and (24), the state weight in eq. (22) can be expanded. Denoting the isospin projection 
coefficients by: 



3=1 

n^-,{^}i 



3=1 



\i,h)(i,h\ 



\i,hw,h\ 



K 



3=1 

fi\ij,{-i? M }) 



(29) 



3=1 



the state weight turns into, by using the factorization of projectors in eq. (3): 



K 



Y\X(P3 



3 = 1 



K Lj 

4 n > } (/, j 3 ) (s„\p PJXn p v \s p ) n n <w>. 

3=1 h=i 



> Hi 



+ C Xc T { p N]} (I, I 3 ) ( Sl \P P jxuPv\s p ) n n % 



K Lj 



(30) 



3=113=1 



,qi,- 



where Kronecker factors S I, and S a 3 ,, , stem from the vanishing of scalar products between single particle states 

with different baryon number and strangeness. Note that, as the parity of particles and antiparticles are the same by 
construction: 



\ s p) = K*«(li)}.{»? w (li)}> = K fc «(li)}.{»7«(li)}> = \ s p) 



(31) 
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therefore, taking into account eq. (31), we rewrite eq. (30) as: 



K 



3=1 



{s L \Ppjxn?v\s P ) 



(32) 



3 = 1 h=l 



The microcanonical weight of a multi-hadronic channel can be then calculated (for spherical clusters) from (32) by 
summing over particle momenta p and polarizations a and averaging over the initial polarization of the cluster. In 
formula: 



I y 



(2J+1) 



K 

nE 



K 



n / d3 p« 3 



(33) 



where if is the number of types, Nj is the number of particles of the type j and the factor rij=i ^l^r 1Sl needed in 
order avoid multiple counting of (anti)symmetric basis tensors when integrating over all particle momenta. 

The integration over kinematical variables, understood in \s p ), gives rise to the same expression as in (15), with 
the difference that now j labeling types and not species. There is also an additional coefficient related to internal 
symmetry (isospin and C-parity): 



n {Jv . } = <o|p y |o>£ 



f[x( Pl ) b 

3 = 1 



4tt 



K.N, 
j=l J n 3 =l J 



N 



P y^ Pn I sni — sin 



t -■- 
2 



K 

n 

3 = 1 



n 



i=l L 



sm[(Sj + 



-\ hiApj)' 



x (i^ } (/,/ 3 ) n n +xxd P N]} (i,h) n n s -^ h ^ 

3 = 1 h = l 3 = 1 h = l 



(34) 



The eq. (34) is the final expression of the microcanonical channel weight for a multi-hadronic channel. This applies 
to completely neutral clusters with C-parity X- F° r charged clusters the microcanonical channel weight can be obtained 
by removing the C-parity projector, i.e. setting \ — in the above equation and multiplying the result by 2. 



IV. NUMERICAL COMPUTATION 



The channel weight in eq. (34) is the basic expression to calculate exclusive channel rates; it cannot be worked 
out analitically, but it can be evaluated numerically. According to eq. (34), the task is indeed twofold: first, the 
computation of the isospin coefficients (see eq. (29) and, second, that of multi-dimensional momentum integrals. The 
sum over permutations can be made with well known algorithms. In order to compute the isospin coefficient we 
have designed a suitable recursive algorithm which is described in Appendix B, while in this section we focus on the 
problem of computing the momentum integrals. 

For properly relativistic particles this problem is known not to have an analytic solution. Previous attempts to 
obtain sufficiently accurate estimates [16] resorted to Monte- Carlo integration and this is the technique we take. For 
the sake of simplicity, we will describe our method for a channel with particles of different species; the generalization 
to identical particles entails a sum over permutations in the integrand function and does not involve special difficulties. 
In this case, for a cluster at rest, the general momentum integral in eq. (34) can be written as: 



J d 3 Pl . . . d 3 PJV( 5 - y s,}j s 3 (y P „^J w({ Pn }) 



(35) 
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where the function W is, up to a constant: 

t>4ir 



W({ Pn }) 







A, • * • 

dtp sin — sm 



N 

n 

n=l 



sin[(5„ + 



sin(|) 



(36) 



N N 

?( s )i„ D-l^/A^ ^ l TTTT.TT T?( s )l 



X II F v(Pn - Ri" 1 WPn) + IIII, [] FviPn + R^WOPn) ■ 



\n=l 



n=l 



In eq. (35) the last momentum variable p^ can be integrated away at once with the 5 3 of momentum conservation, 
yielding: 

. / N-2 \ 

J d 3 Pl . . . d 3 PA r_ 1( 5 f M - £ ™ - \Jp%-i+ m N-i - ^\-Pn-2-Pn-i\ 2 +m 2 N J W({p n }) 



where: 



N-2 



N-2 



E P" • 



n=l 



Note that in eq. (37), for simplicity, we have used the same symbol for W({p n }) keeping in mind that indeed, after 
the integration, this is a different function as it only depends on the set of momenta pi, . . . , pjv-i- Denoting with 
p^! and p^f2i the two zeroes of the argument of the S in eq. (37), this can be rewritten as: 



d 3 pi...d 3 PA r_i 



£n£n-i 



Pn-i(£n-i + £a0 + £n-iPn-2 ■ Pn-i 

6(PN-1 - P ( N-i) + 6(PN-1 ~ P ( N-1)} W({ Pn }) 



(37) 



where Pn-i is the versor of pjv_i and p^2i an d p)^\ fulfill the following equation: 



(B) 



N-2 



(pN- 



+ - 



N-2 



n=l 



p2 

N—2 



/ N-2 



Pat_i + (pjv-i • Pjv-2) 



N-2 



p2 

r N-2 



+ m 2 N - m 2 N _ 1 



■ 2 2 
+ m N — m 



N-l 



71=1 



N-2 



'N- 



-1 M -E 



n=l 



. 



PJV-1 



n=l 



If wc now let: 



-{A) (A) 
L N t N-l 



Pjv-i ( e iv-i + e iv ') + 4f-i p iv-2 • Pjv-i 



(38) 



(39) 



where: 



£ (A) - 
N — 



m 



N-l 



Pw-iPw-i + Pn-2 +m 2 N and = y (p^-i) 

and similarly for p N _n the eq. (37) can be finally written as: 

J d 3 Pl . . . pjv-2 dOjV-! [$ A W({p„} A ) + $ B ^({p„} B )] (40) 

where fijv-i is the solid angle of pjv-i and VF({p„} A ) stands for the function W in eq. (36) evaluated on the momenta: 

f ^A f (A) , ( (A) - . ry \ \ 

\Pni = ( Pi, • • • , PiV-2,Pjv_lPiV-l, — ( P N _iPN-l + ^ N -2 J J ] 

and similarly for W({p n } B ). 

In order to calculate eq. (40) we have to perform a 3 AT — 4-dimensional momentum integration plus one further 
integration over the angle ip hidden in the function W. Overall, this is a 3iV — 3 dimensional integration which is 
carried out by using an importance sampling Monte-Carlo method, designed to achieve the best performance. 
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A. Importance Sampling 



The importance sampling method is a well-known method to perform Monte-Carlo integration. The idea is to find 
an auxiliary function g(x) which is at the same time easy to sample and as similar as possible to the integrand function 
to maximize efficiency. Each random extraction in x is then weighted by the ratio f{x)/g(x) and an estimator of the 
integral is given, after N extractions, by 

JV 

^2f{x l )/g{x l ) 

i=l 

In our case, we have to extract 3iV — 3 variables: N — 2 momenta, N — 1 solid angles, and one further angle ip. 
We extract all angles according from a flat distribution while for momenta, our method is based on the use of the 
asymptotic limit (N — > oo) of the integrand; this is a known one, as for N — > oo the microcanonical iV-body phase 
space should converge to its canonical ensemble limit, which consists - in the Boltzmann limit - of a factor exp[— s/T] 
for each particle. Therefore, for each particle, we expect: 



p 2 dpc- £/T dt^t(t + 2m i )(t + nij)c- t/T , (41) 

for its kinetic energy distribution. The temperature T in eq. (41) is a parameter to be chosen to minimize the difference 
between integrand and auxiliary function. We have calculated it, along with a set of chemical potentials associated 
to each conserved charge (electric, baryonic and strange) by equating the total energy, momentum and charges in the 
microcanonical ensemble with their average value in the grand-canonical ensemble, which is precisely the saddle-point 
equation governing the asymptotic expansion of the microcanonical partition function [14]: 

M = T 2 A^. (TK <WT ( 42 ) 

3 

3 

where M is the cluster mass and is the set of abelian charges of the species j; /i is a set of chemical potentials and: 

Zj (T) = (2 Sj + l)^ / d 3 pexp(-e,/T) = t&^mjTK* (^) (43) 

where K2 (nij/T) is the modified Besscl function of the 2 nd kind. 

Unfortunately, the function (41) is not a practical auxiliary function as it cannot be sampled fast enough. Its integral 
cannot be inverted analytically and there are not optimized sampling algorithms either in the ultra-relativistic or non- 
relativistic limit where it can be approximated as: 

t » rrij -> - A(f e" t/T ) (44) 
t < m 3 - A(Vt c~ t/T ) 

where A is a normalization constant. We have then replaced the function (41) with the mathematical j3 function: 

r(a)r(6) 



which can be very efficiently sampled [17] with a fast rejection algorithm [18] . In our case, x is set to be the ratio t/t 



where t is the kinetic energy of each particle and t max is the highest available kinetic energy, that is M — m i' 

p (*).^(*)~'( 1 -<)'- iT ,46, 



max • 



tmax/ T(a)T(b) \t 

The constant b is chosen to be b = 1 + t max /T where T is obtained from solving eqs. (42). This choice is dictated by 
the fact that in the limit t max /T » 1 this distribution is very close to an exponential, like in eq. (41): 

l im (l-^) ^-e-*/' 
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The only problem of the distribution (46) is that it does not match properly the actual one when t is large, i.e. when 
t ~ i m a X ; indeed, eq. (46) goes to zero too rapidly and this makes the algorithm less effective in sampling the region 
in multi-dimensional momentum space where one particle in the channel carries most of the available kinetic energy. 
However, such contributions are relevant only for channels with few light particles and one much more massive. 
As has been shown in [17], a good choice for the parameter a is: 

a-l = \ + \ C - 2 ^ (47) 

where m,j is in GeV. The empirical dependence of a on particle mass is such that a takes the values 3 and 1.5 
respectively in the ultrarelativistic and non-relativistic limit, according to (44). 
Summarizing, the steps of the Monte-Carlo integration algorithm are as follows: 

• extraction of the angle ip with a flat distribution; 

• extraction of N — 1 solid angles with a flat distribution; 

• extraction of N — 2 kinetic energies according to (3 distribution; 

• evaluation of the modulus Pn-i solving eq. (38); 

• evaluation of = — ^2n=i Pn- 

We note that eq. (38) can have either one real solution or two real solutions or none. In case of two solutions 
the integrand is evaluated on both and averaged, while for no solutions the whole extraction procedure is repeated. 
We also take into account, for channels with resonances, their mass broadening. In fact, masses are also randomly 
extracted, at each step of the Monte-Carlo integration, according to a relativistic Breit-Wigner distribution: 

B r (m) = -!- 7 — —— (48) 

y ' 2tt (m-m r ) 2 +T2/4 v ; 

where m r is the central mass value and T r the width. 



B. Resonance decay channels 

As has been discussed at the end of Sect. II, for each exclusive channel we assume the hadron-resonance gas picture, 
which prescribes that the probability of an exclusive channel with some set of final hadrons {Nj} is the sum, with 
suitable weights, of the probability of all possible channels with hadrons and resonances decaying - in a single or 
multiple steps - into the set {Nj}. The weights are given by products of branching ratios of parent hadrons and 
resonances. Therefore, given a final channel {Nj}, the first task is to find all channels which may have generated it, 
knowing the decay modes of all hadrons and resonances. The search of all parent channels is a multi-step recursive 
problem, in that many generation steps can occur. If we denote by {Nj}^) a channel which can directly decay into 
the channel {AT,}, by {Nj}^ 2 ) a channel which can directly decay in {Nj}^ and so on, one has to find all possible 
decay trees like those shown in fig. 2. In view of the large number of resonances, this task is a hard combinatorial 
problem and a suitable recursive algorithm has been designed to solve it. 

Our method is to list all subsets of particles in the channel { Nj} and check whether there is one or more resonances 
decaying in each subset. For instance, for a channel with four particles labelled with 1, 2, 3, 4, possible subsets are: 



(1)(2,3,4) 


(2,3)(1)(4) 


(2)(1,3,4) 


(2,4)(1)(3) 


(3)(1,2,4) 


(3,4)(1)(2) 


(4)(1,2,3) 


(1,2)(3,4) 


(M)(3)(4) 


(1,3)(2,4) 


(1,3)(2)(4) 


(M)(2,3) 


(M)(2)(3) 


(1,2,3,4) 



Each resonance decaying into a subset gives rise to a possible upper-level parent channel made of the resonance 
replacing the given subset and the remaining hadrons of the original channel. For each parent channel found, the 
procedure is iterated. At the uppermost level, possibly one has a channel with only one massive resonance having the 
same quantum numbers as the initial state. 
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FIG. 2. Examples of possible decay trees for a four particles channel. Circles encompass decay products of the particle at 
higher level. 




FIG. 3. Examples of possible equivalent decay trees for a four particles channel. 



The complexity of this combinatorial problem is considerable. Actually, the number B n of subsets of a set of n 
integers is known as the Bell number which is given in terms of a recurrence relation: 

S„ +1 =X>fcQ- (49) 



This number grows quickly as a function of n, so that listing all possible parent channels becomes practically impossible 
for channels with more than ~ 8 particles (B& = 4140 and Bg = 21147), also in view of the very large number of 
resonances in the hadron spectrum. 

Indeed, this recursive search algorithm can give rise to multiple counting. This happens if there are identical 
particles in the channel simply because some subsets are actually the same. For instance, if there are four particles, 
the subsets (1,2)(3,4) and (1,3)(2,4) are clearly equivalent if particles 2 and 3 are identical. This kind of double 
counting is quite easy to get round algorithmically a priori. Yet, double counting may occur even if particles in the 
channel are all different. Consider for instance two decay trees like those shown in fig. 3 with, on the left, the subset 
(a)(b)(cd) and on the right the subset (ab)(cd). In fact, both configurations can stem from the same channel e, / at a 
different level, so this parent channel may appear twice in our parent channel search. This kind of multiple counting 
is avoided by re-checking a posteriori the full list of parent channels found. 

Finally, the probability P{Nj} of observing a final channel {Nj} can be written as a finite sum over all parent 
channels: 

P{Nj} °C ^{Nj} = ^{N 3 } + BR (l)^{JV 3 } (1) + BR (2) BR (l)^{A' J }(2) + • • • 

+ BR {iy n {Nj}(i)l + BR (2) ,BR (1) ,n W}(2)( + . . . (50) 

where BR(j) is the product of branching ratios of particles in the channel {Nj}u\ decaying into particles in the channel 
{-/Vj-}(j_i) and where ^{Nj} is, by definition, the total channel weight, including contributions of parent channels. 
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V. ANALYSIS OF e+e COLLISIONS AT LOW ENERGY 

As has been mentioned in the Introduction, the rates of exclusive hadronic channels can be measured only in low 
energy collisions (say < 5 GeV) because the large multiplicity of the final state at high energy makes a full identification 
of particles impossible. There have been in the past some attempts to reproduce hadron multiplicities and some 
multi-pion(kaon) differential cross sections in low energy e + e~ collisions [2] by using statistical-thermodinamical or 
statistical-inspired models. Yet, in none of those calculations the full set of conservation laws has been taken into 
account, because of the lack of a proper formulation of the microcanonical ensemble with intrinsic angular momentum 
and the involved numerical calculations. As we will show in Subsect. VA, this is a serious drawback because, when 
dealing with exclusive channels, all conservation laws play a major role. In fact, we now have all needed ingredients 
to make a proper test of the statistical hadronization model with exclusive channels: the statistical weight of multi- 
particle channels including exact energy-momentum and intrinsic angular momentum conservation [12] (formula (34) 
and a sufficient computing power. 

As discussed in Sect. II, at low energy the formation of a single cluster at rest in the centre-of-mass frame of an 
c + c~ collision is assumed. Its mass will therefore coincide with \fs and the other quantum numbers will be those of the 
initial state. Particularly, in e + e~ collision, the hadron production is dominated by the diagram with an intermediate 
virtual photon (see fig. 1), so that the hadronizing cluster is assigned with a spin, parity and C-parity J PC = 1 . 
On the other hand, isospin is not conserved in electromagnetic interaction and it is therefore unknown; in the Vector 
Dominance Model (VDM) this depends on the coupling of the photon to different resonances, but we will be working 
in a mass region above 2 GeV, far from known resonance region (see discussion in the following). Therefore, we will 
assume an unknown statistical mixture of I = and 1=1 initial state, neglecting interference terms, and introducing 
a free parameter I such that for the mixed state: 

Io|0,0)(0,0| + (1-I )|l,0)(l,0|, 

Finally, the geometry of the cluster needs to be fixed. We assume a spherical shape and a volume given by: 

v=™ = ^ (51) 
p p 

where M is the mass and p the energy density; this is taken to be a free parameter to be determined by comparing 
the model with the data. 

Motivated by observations concerning hadron abundances at high energy, we allow deviation from the full statistical 
equilibrium of channels involving particles with strange valence quarks. This is done by introducing in the analysis 
the extra-strangeness suppression parameter 75 [19]. For its definition here to be in agreement with the formulae 
of inclusive multiplicities of hadrons in the canonical and grand-canonical ensembles, one just needs to multiply the 
microcanonical weight of a channel by j s 3 , Sj being the number of valence strange quarks of each particle: 

n {Nj} (52) 

The 7s factor also applies to neutral mesons with hidden strange quark content like 77, <j> etc. Since the wavefunction of 
such particles is in general a superposition like C u uu + Cddd-|-C s ss with |C U | 2 + |C d | 2 + |C S | 2 = 1, only the component 
ss of the wavefunction is suppressed, i.e. we multiply by: 

[|C S | 2 7§ + (1-|C S | 2 )] . 

for each neutral meson. To calculate |C S | 2 , we have used mixing angles quoted by the Particle Data Book [20]. 

The branching ratios, masses and widths of hadrons and resonances needed to calculate the exclusive channel rate 
according to formula (50), have also been taken from the latest issue of the Particle Data Book [20]. All hadrons up to 
a mass of 1.8 GeV for mesons and 1.9 for baryons have been included for the generation of parent channels. It is now 
appropriate to briefly discuss the possible contribution of single resonance decay (off-peak) to the hadron production 
in c+c - collisions at low energy. In terms of the diagrammatic description in fig. (2), these contributions correspond 
to the highest ancestor of the channel in the decay tree, being a single resonance with the same quantum numbers 
as the initial state. This contribution, if relevant, cannot be subtracted away from the experimental data given the 
poor knowledge of resonances above 1.8-1.9 GeV mass. If one assumes that resonances can be identified with clusters 
decaying statistically, then this contribution should be somehow taken into account within the SHM calculation itself. 
On the other hand, if we refrain from this assumption, we must move sufficiently far from the energy region where 



n 



K 

n 



,7s 
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(GeV) 


p (GoV/fm 3 ) 


7s 


C* P7S (GeV/fm 3 ) 


Io 


A (nb GoV 4 ) 


X 2 /dof 


2.1 


0.24 ±0.17 


0.66 ±0.22 


0.011 


0.17 ±0.03 


0.028 ± 0.002 


93.4/16 


2.2 


0.36 ± 0.20 


0.86 ±0.22 


0.023 


0.15 ±0.03 


0.036 ± 0.002 


82.6/14 


2.4 


0.44 ±0.30 


0.78 ±0.36 


0.024 


0.23 ±0.04 


0.017 ±0.001 


55.4/17 


2.6 


0.56 ± 0.36 


0.62 ±0.47 


-0.009 


0.53 ±0.07 


0.017 ±0.002 


44.9/12 



TABLE I. Summary of the fit results to multi-hadronic exclusive channels at different centre-of-mass energies. Also shown the 
correlation coefficient of p and 73. 




FIG. 4. Upper panel: best fit energy density p as a function of centre-of-mass energy. Lower panel: best fit 7s values. Lines 
are drawn to guide the eye. 



narrow resonances lie in order to minimize their impact on the cross section. Furthermore, we do not want to get 
over the charm quark production threshold and this constrains the energy interval to about 2-3 GeV. 

Much data in this energy interval has been lately provided by the BABAR experiment which has measured the 
cross-sections of several multi-hadronic channels in c + c~ collisions at several centre-of-mass energies with the method 
of initial state radiation. We have chosen four energy points, that is y/s = 2.1, 2.2, 2.4 and 2.6 GeV and added to 
the available BABAR measurements older measurements performed by experiments at e + e~ colliders run at the same 
centre-of-mass energies and collected in a nice review paper [21]. In fact, most centre-of-mass energies were near, but 
not exactly, those values. Therefore, for each chosen energy, we have interpolated cross sections by making a simple 
mean of cross sections and errors measured at adjacent energy values. For each energy point, the thereby obtained 
cross-section values and errors from different experiments have been further averaged according to the weigthed average 
method used by the Particle Data Group [20], including error rescaling by x 2 /dof in case of large discrepancy. The 
full set of channels can be read through in tables VI, VI, VI, VI in Appendix C. 

In order to compare the calculation with the data of exclusive channels rate, given in terms of a cross section, we 
have introduced a normalization free parameter A(y/s): 



<7 {N j } = a (Vs)u{n 3 } 



(53) 
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Finally, we have fitted all available measurements of exclusive channels rates at a given energy to the SHM with four 
free parameters: p, 75, A and Iq. The fit minimizes the x 2 : 

r . l ^oxp "r ^thoo 

L J J measured 

where the sum runs over measured channels; A exp is the experimental error and A t hoo is the theoretical uncertainty on 
the cross sections respectively. The latter is the sum in quadrature of the statistical error, owing to the finite statistics 
in Monte-Carlo integration and the systematic error stemming from the uncertainty on branching ratios of resonances. 
This error has been estimated at each step of the x 2 minimization by varying the branching ratios by their errors quoted 
in the Particle Data Book (or by making an educated guess if missing) and calculating the difference between the 
theoretical value of the probability (50) before and after the variation. This additional uncertainty generally overcomes 
the experimental error for ligth particle channels with > 2 particles, as the number of contributing channels with 
resonances is large (as it can be seen in tables VI, VI, VI and VI in Appendix C). 

The fit procedure is as follows: for each energy point a two-dimensional grid in the parameters p and 75 is set, 
with 50 divisions in each direction and range [0.04-2] GeV/fm 3 , [0.02-1] respectively. At each point of the grid, a 
minimization of the x 2 is carried out to determine the parameters 1$ and A. The point grid where the minimum 
among all minima lies has been taken as the best fit. The error on the parameters p and 75 has been estimated 
graphically from the contour x 2 = Xmin + 1 ( see n S- 5)- 




FIG. 5. Contour plot of the X 2 for yfs = 2.4 GeV. 



The fit results are shown in summarized in table V and in fig. 4. The comparison of the fitted rates with the data 
is shown in figs. 8, 9, 10, 11 and in tables VI, VI, VI, VI in Appendix C for y/s = 2.1, 2.2, 2.4, 2.6 GeV respectively. 



A. Test of conservation laws 



To appraise the role of the various constraints and the features of the statistical model, it is worth making the fit 
by switching off the conservation laws in turn. We have therefore minimized the x 2 with the same aforementioned 
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FIG. 6. Fits to exclusive channels in e + e~ collisions at 2.4 with conservation laws turned off. The data points, shown as full 
dots, are weighted averages of available measurements (see text) while the fitted points are shown as open squares. Panel (a): 
fit with only energy- momentum conservation. Panel (b): fit without angular momentum conservation. 
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procedure, at a single energy point y/s — 2.4 GcV, in two different modes. 

In the first mode we have kept only the conservation of energy and momentum, disregarding angular momentum, 
parity and internal symmetries. The result is shown in fig. 6(a), where a consistent worsening of the fit quality can 
be seen (x 2 /dof = 133/17), especially looking at the residual distribution. The fit quality significantly improves 
by turning on the internal symmetries (% 2 /dof = 74.5/17). Finally, restoring the angular momentum and parity 
conservation, one obtains the best fit shown in fig. 10 (\ 2 /doi = 55.4/17). This further improvement indicates that 
angular momentum conservation plays an important role and this was indeed expected as this is a very effective 
mechanism in modulating the rate of a channel for a spacially extended source (the well known centrifugal barrier 
effect), which is one of the crucial assumption of the statistical model. 
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FIG. 7. Normalized microcanonical weights of the channel n it as a function of the cluster spin J in the Boltzmann statistics 
(upper panel) and quantum statistics (lower panel) case. The cluster has been taken spherical in shape with radius R — 4 GeV -1 
and a mass M = 5 GeV. Isospin, parity, C-parity conservation are turned off. 



To highlight this effect, we have calculated the microcanonical channel weight for the channel 7r°7r° for different spin 
J of the cluster. Assuming a spherical shape and switching off parity and internal symmetries (isospin and charge 
conjugation) we have obtained the normalized microcanonical channel weights shown in fig. (7), for a cluster mass 
M of 5 GeV and a radius of R = 4 GeV -1 . As expected, the maximal value of the microcanonical channel weight 
is located around the angular momentum J <~ pi? ^ 10, p being the momentum of the pion. This calculation also 
served as a consistency check for our numerical code as the sum of all microcanonical channel weights for all J was 
found to reproduce the simple two-body phase space with only energy-momentum conservation. 



VI. DISCUSSION AND CONCLUSIONS 



Although the fit quality is not perfect in terms of statistical test (see table V), we can fairly conclude, looking at 
figs. 8, 9, 10, 11 that the statistical hadronization model is able to satisfactorilty reproduce most exclusive multi- 
hadronic channels measured in e + c~ collisions at low energy. Especially at 2.4 GeV, all measured channel rates lie 
within 2.5 standard deviations from the model values, which is quite remarkable taking into account the obvious fact 
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FIG. 8. Upper panel: exclusive cross sections of various channel in e + e~ collisions at 2.1 GeV. The data points, shown as full 
dots, are weighted averages of available measurements (see text) while the fitted points are shown as open squares. The lines 
have been drawn to guide the eye. Lower panel: fit residual distribution. 



that exclusive channels are a very stringent test for any model, certainly much more than inclusive multiplicities, and 
that the fits were done with only 4 free parameters. 

To fairly judge the quality of the agreement between model and data, it is worth keeping in mind that the analysis 
we have presented in this work still relies on several approximations, so that one may hope that a more thorough 
calculation will result in a better agreement between model and data. The main approximation introduced in this 
calculation is the assumption of validity of the hadron- resonance gas model for exclusive channels and at finite volume, 
what is in fact granted only for fully inclusive quantities and in the thermodynamic limit. Moreover, the assumed 
symmetry SU(2)j so £g)U(l) s trange may not be a sufficiently realistic scheme. 

The statistical model nicely matches one of the main features of the data, namely the gross dependence of the 
channel rate on the number of particles TV; this kind of hyerarchy can be clearly observed in all fits (see e.g. fig. 10) 
and the SHM is able to reproduce this behaviour because of the approximate dependence of the rate on V N , as in 
formula (18). Furthermore, the more conservation laws are included, the more the fit improves and the more the 
model predictions approach the actual rates. This is a certainly a good point for the SHM, especially as far as angular 
momentum is concerned, as has been discussed in Subscct. VA. 

Overall, the most interesting outcome of the analysis are the values of the fitted energy density p and strangeness 
suppression parameter 7s, shown in fig. 4, around 0.5 GeV/fm 3 and 0.7 respectively. These values are essentially 
the same obtained with the analysis of inclusive hadronic multiplicities at high energy [3] 4 and this confirms the 
consistency of the statistical approach to hadronization. While the origin of extra-strangeness suppression is not clear 
(an interesting proposal was put forward in ref. [22]) the idea of hadronization as a process occurring at a critical 
energy density which uniformly populates the available phase space is certainly reinforced by the observation that 
this seems to happen at universal values of the parameters, at high as well as at low energy. 



4 For a neutral hadron-resonance gas in the canonical ensemble, the energy density of 0.4 GeV/fm 3 approximately corresponds to a 
temperature of 160 MeV 
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FIG. 9. Upper panel: exclusive cross sections of various channel in e + e~ collisions at 2.2 GeV. The data points, shown as full 
dots, are weighted averages of available measurements (see text) while the fitted points are shown as open squares. The lines 
have been drawn to guide the eye. Lower panel: fit residual distribution. 



The interpretation of the statistical equilibrium in hadronization is an open issue and there are several proposals. 
One of the authors (F.B.) favours the idea of a quantum-chaotic effect (known as Berry's conjecture) owing the the 
strong non- linearity of QCD in the non-perturbative regime. 
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APPENDIX A 



The charge- conjugation operator C, when acting on a charged light-flavoured meson belonging to an isotriplet may 
generate an arbitrary phase. Fixing this phase is essential for our definition of the projector (23). For this purpose, 
let us define the operator G as: 

G = C e i7r?2 (55) 

which is known as G-parity operator. When applied to a state with third component of isospin I3, this amounts to 
first flipping J 3 -4- — -Z3 and reversing the process, — ^3 -> I3. Therefore, an eigenvector of I3 is also an eigenvector of 
the G-parity operator. For an isotriplet like pions': 

G|tt ± > = — Cj^) (56) 
G\n°) = ±k°) . 
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because e l7rl2 results in a factor —1 when applied to any pion state. Now, since \n°) is an eigenstate of C with 
eigenvalue +1, we define: 

CIO = -Itt*) . (57) 
so as to G-parity to yield the same eigenvalue for all members of the isotriplet. 



APPENDIX B 



One of the important steps in the numerical evaluation of the microcanonical channel weight is the computation of: 



r{N 3 } 



(I, h) 



K 



\i,h)(i,h\ 



K 



(58) 



(59) 



We describe here a method to compute the more general expression: 

{(i\m^i!\---(i\i^\}\ij 3 )(ij 3 \{\i\ii')\i^ii')---\i N ,ir)} 

where 

\I\II)\IMI)---\I N ,I 3 N ) (60) 

is a generic multi-particle isospin state. A closed analytical formula for (58), as a finite sum, has been found in 
ref. [23]. We have not used that formula, yet, in the following, we will closely follow the notations therein. 

In order to calculate the projection of the state (60) onto the subspace with total isospin \II 3 ), one should use a 
system of base- vectors in isospin space where / is diagonal. The choice of such a basis is equivalent to the choice of a 
coupling scheme for the Ii [23]; we can choose, for instance, a base where: 

(I 12 ) 2 = (I 1 + I 2 ) 2 (61) 

( 7 123 ) 2 = (7 l +7 2 + / 3 ) 2 



are diagonal. 

A base vector for this scheme is denoted by: 



I/ 12 ,/ 123 , ...,I,h 



(62) 



where 1,1, etc., are known. Both sets (62) and (60) are complete and they are connected by a unitary transformation: 
|/ 12 ,/ 123 ,...,/,/ 3 ) (63) 

= £ {(/M3 1 |...(^,^|}|/ 12 ,/ 123 ,...,/,/ 3 ){|/ 1 ,/ 3 1 >---|^/3 A '>} 

The coefficient { (J 1 , I 3 \ ■ ■ ■ (I N , I 3 \} \I 12 , 1 123 , . . . ,I,I 3 ) is called a recoupling coefficient and is a product of Clebsch- 
Gordan coefficients: 



Sit 1 t 1 \ it n t n \\\t 12 t 123 t t \ - r ,/1/2/12 r i 12 n - 1 i n i 



By using eq. (64) one can rewrite (59) as: 



{{l\ll\...{I N ,I»\}\I,h){I,h\{\l\ll>)...\I N ,I»>)} 

J2 {{i\n\---{i N ,i 3 N \}\i' 2 ,i^,...,i,h) 

/12/123 J12...N — 1 

x(l' 2 ,I^,...,I,I 3 \{\l\l^')---\I N ,li")} 



(64) 
(65) 



Wi/f (/i+/f ) ' ' ' (/i+/f+.../ 3 ~- 1 )/« 

J12J123 J12...N — 1 



^z 1 / 2 / 12 n i 12 n - 1 i n i 

xo /l >q '(ii ') ' ' ' <+/f >)i? >i 3 
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of course, if one of the Clebsch-Gordan coefficient is vanishing, the corresponding term in the previous sum is vanishing 
too. 

We have implemented a recursive numerical method to calculate the above expression. By iterating the coupling 
scheme in (61) one can build a tree diagram, shown in fig. (12) from left to right, where for each recoupling step, the 
highest value of the resulting isospin is put on top and the other values are sorted in decreasing order. 




FIG. 12. Tree diagram for the calculation of isospin matrix elements. 



To each branch of the tree a numerical coefficient is associated, which is the recursive product of the Clebsch-Gordan 
coefficient corresponding to the decomposition which generated that branch and the numerical coefficient of the gener- 
ating branch. At the righmost end of the tree, one finds all possible values of the global isospin which can be obtained 
by coupling the set I 1 , . . . ; I N , 1$ and the corresponding recoupling coefficients (I, hll 1 , 13 , 1 2 , If , . . . I N , 

The main advantage of this algorithm resides in the possibility of a simple recursive implementation. Moreover, the 
CPU time needed for the tree exploration can be easily reduced by simply switching off the recursion for branches 
with vanishing Clebsch-Gordan coefficients. 



APPENDIX C 



In this Appendix we have collected the tables with best fit values and the experimental values of the cross sections 
of exclusive channels. 
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v^= 2.1 GeV 


channel 


itsm (nb) 


crexp. (nb) 


References 


7T + 7r _ 


0.229 ± 0.003 ± 


0.177 ±0.036 


[21, 24] 




0.133 ±0.002 ±0.03 


0.385 ±0.11 


[21, 25] 




5.93 ±0.06 ±0.5 


6.2 ±2.57 


[21] 




3.96 ±0.02 ±0.4 


4.18 ±0.45* 


[21, 26] 


_+^-_o_o_o 

7T 7T 7T 7T 7T 


0.809 ±0.02 ±0.05 


4.03 ±2.41 


[21] 


7r + 7r _ 7T + 7r _ 7r 


1.46 ±0.01 ±0.07 


1.46 ±0.4* 


[21, 27] 


?77r + 7r _ 


0.758 ±0.02 ±0.09 


0.36 ±0.11 


[21, 27] 


W7T 


0.261 ± 0.004 ± 


0.68 ±0.19 


[21] 




0.0673 ± 0.0007 ± 


0.18 ±0.11 


[21] 


0J-K + 1T~ 


0.478 ± 0.002 ± 0.001 


0.3 ±0.09 


[21, 27] 


riir + n~ iv + n~ 


2.18 ±0.02 ±0.7 


0.915 ±0.2 


[27] 


K+K-tv 


0.148 ±0.005 ± 0.0009 


0.12 ±0.042 


[21, 28] 


K+K-7V°7T° 


0.242 ± 0.002 ± 0.006 


0.9 ±0.134 


[28] 


K+K-TV+TY- 


0.929 ±0.003 ±0.05 


3.05 ±0.5* 


[21, 28] 


K+K-TY+TY-n 


0.6 ±0.006 ±0.03 


0.575 ±0.13 


[27] 


n<t> 


0.0206 ± 0.0002 ± 


0.425 ±0.17 


[29] 


w/o 


0.151 ±0.0007 ±0 


0.39 ±0.11 


[27] 




0.479 ± 0.002 ± 


0.63 ±0.076 


[21] 


nn 


0.48 ± 0.002 ± 


1±0.28 


[21] 




0.0331 ± 6e - 05 ± 


0.019 ±0.007 


[28] 


Not included in the fit 


r]'ir + -K 


0.371 ± 0.005 ± 0.008 


0.17 ±0.07** 


[27] 




0.127 ± 0.0009 ±0 


0.395 ±0.065** 


[28] 




0.324 ±0.002 ±0.01 


1.07 ±0.08** 


[28] 


/l7r + 7T~ 


1.06 ±0.005 ±0.2 


0.79 ±0.24** 


[27] 




0.24 ± 0.001 ±0 


0.49 ±0.09** 


[27] 


4>fo 


0.089 ± 0.0002 ±0 


0.41 ±0.08** 


[28] 



* Errors have been rescaled because of discrepancies between different experiments. 
** Errors are statistical only. 



TABLE II. Comparison between fitted and experimental cross sections of exclusive multi-hadronic channels in e + e~ collisions 
at y/s — 2.1 GeV. The errors on theoretical cross sections oshm are the Monte-Carlo integration statistical error and the 
error owing to uncertainty on branching ratios of resonances contributing to the channel, respectively. The experimental cross 
section values have been obtained by averaging available measurements, see text. 
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y/s = 2.2 GeV 


channel 


(tshm (nb) 


&ex P . (nb) 


References 




0.188 ± 0.002 ±0 


0.101 ±0.0359 


[24] 




0.0971 ±0.003 ±0.02 


0.395 ±0.17 


[25] 




4.53 ±0.03 ±0.4 


4.2 ± 1.2 


[21] 




3.34 ±0.1 ±0.3 


3.73 ±0.2 


[21, 26] 


7r + 7r _ 7r + 7r _ 7r° 


1.21 ±0.01 ±0.06 


1.25 ±0.26 


[27] 


?77r + 7r _ 


0.712 ± 0.007 ±0.1 


0.48 ±0.12 


[27] 




0.207 ± 0.003 ±0 


0.41 ±0.2 


[21] 




0.355 ± 0.002 ± 0.0008 


0.2 ±0.09 


[21, 27] 


?77r + 7r _ '7r + 7r _ 


2.57 ± 0.03 ± 1 


1.09 ±0.215 


[27] 


K+K-tv 


0.166 ±0.02 ±0.001 


0.146 ±0.05 


[29] 


K+K-7v ir 


0.34 ± 0.002 ± 0.008 


0.57 ±0.1 


[28] 




1.29 ±0.005 ±0.06 


2.07 ±0.16 


[21, 28] 


K+K-TY+n-Tv" 


1.34 ±0.008 ±0.07 


0.54 ±0.13 


[27] 


V<f> 


0.0285 ± 0.0003 ± 


0.235 ± 0.126 


[29] 


w/o 


0.0976 ± 0.0004 ± 


0.17 ±0.091 


[27] 




0.0134 ± 0.0003 ±0 


0.02 ±0.02 


[27] 




0.34 ±0.001 ±0 


0.454 ± 0.052 


[21] 




0.0648 ± 0.0003 ± 


0.0635 ±0.015 


[28] 


Not included in the fit 




0.341 ± 0.005 ± 0.007 


0.101 ± 0.052* 


[27] 


</>7T~ + 7T~ 


0.139 ± 0.001 ±0 


0.27 ±0.055* 


[28] 


/r *(892)if-7r+ 


0.508 ± 0.003 ± 0.02 


0.62 ± 0.062* 


[28] 


/l7T + 7T~ 


1.02 ±0.005 ±0.17 


0.915 ±0.24* 


[27] 




0.3 ± 0.002 ±0 


0.39 ±0.08* 


[27] 


0/0 


0.087 ± 0.0003 ± 


0.295 ± 0.065* 


[28] 



* Errors are statistical only. 



TABLE III. Comparison between fitted and experimental cross sections of exclusive multi-hadronic channels in e + e~ collisions 
at yfs = 2.2 GeV. The errors on theoretical cross sections oshm are the Monte-Carlo integration statistical error and the 
error owing to uncertainty on branching ratios of resonances contributing to the channel, respectively. The experimental cross 
section values have been obtained by averaging available measurements, see text. 
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= 2.4 GeV 


channel 


<?SifM (nb) 


&exp. (nb) 


References 


TV + TY~ 


0.0716 ± 0.0009 ± 


0.105 ±0.0323 


[24] 




0.0611 ±0.001 ±0.01 


0.26 ±0.13 


[25] 




2.25 ±0.03 ±0.4 


2.7 ± 1 


[21] 


7T + 7r~7T + 7r _ 


1.68 ±0.02 ±0.3 


2.12 ±0.2 


[26] 


7r + 7r _ 7r + 7r _ 7r° 


1.1 ±0.009 ±0.1 


0.975 ±0.23 


[27] 




0.289 ±0.01 ±0.06 


0.145 ±0.06 


[21, 27] 


W7T 


0.0785 ± 0.001 ± 


0.38 ±0.2 


[21] 




0.294 ± 0.002 ± 0.04 


0.09 ±0.09 


[21, 27] 


?77r + 7r _ '7r + 7r _ 


1.64 ±0.02 ±0.6 


0.865 ±0.19 


[27] 




0.0739 ± 0.001 ± 0.0009 


0.111 ±0.036 


[29] 




0.188 ±0.002 ±0.007 


0.335 ± 0.07 


[28] 


K+K--K+-K- 


0.912 ±0.004 ±0.1 


1.3 ±0.167* 


[21, 28] 


K+K-TY+TY-n 


1.39 ±0.006 ±0.2 


0.785 ±0.15 


[27] 


n<t> 


0.0127 ± 0.0002 ±0 


0.06 ± 0.0534 


[29] 


w/o 


0.0558 ± 0.0003 ± 


0.06 ±0.05 


[27] 




0.0604 ±0.0007 ±0.01 


0.025 ±0.05 


[27] 




0.129 ± 0.0005 ±0 


0.146 ±0.03 


[21] 


nn 


0.128 ± 0.0005 ±0 


0.69 ±0.29 


[21] 


K+R-R+K- 


0.0374 ± 0.0001 ± 0.0002 


0.077 ±0.018 


[28] 


AA 


0.0464 ± 0.0001 ±0 


0.133 ±0.033 


[30] 


AS 


0.0883 ± 0.0003 ± 


0.04 ±0.021 


[30] 


Not included in the fit 


?/7r + 7r 


0.125 ±0.003 ±0.003 


0.045 ± 0.04** 


[27] 




0.112 ± 0.0006 ±0 


0.125 ±0.04** 


[28] 


^(892)JsT-7r + 


0.245 ± 0.002 ± 0.008 


0.425 ± 0.045** 


[28] 


/l7r + 7T _ 


0.462 ± 0.002 ± 0.07 


0.557 ±0.17** 


[27] 




0.229 ± 0.001 ± 0.044 


0.23 ±0.06** 


[27] 


0/0 


0.0395 ± 0.0002 ± 


0.135 ±0.045** 


[28] 



* Errors have been rescaled because of discrepancies between different experiments. 
** Errors are statistical only. 



TABLE IV. Comparison between fitted and experimental cross sections of exclusive multi-hadronic channels in e + e _ collisions 
at \fs = 2.4 GeV. The errors on theoretical cross sections oshm are the Monte-Carlo integration statistical error and the 
error owing to uncertainty on branching ratios of resonances contributing to the channel, respectively. The experimental cross 
section values have been obtained by averaging available measurements, see text. 
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y/s = 2.6 GeV 


channel 


fSHM (nb) 


crexp. (nb) 


References 


n + n~ 


0.0354 ± 0.0005 ± 


0.0359 ±0.0162 


[24] 




0.114 ± 0.002 ±0.03 


0.175 ±0.105 


[25] 




0.965 ± 0.005 ± 0.2 


1.44 ±0.16 


[26] 




3.52 ±0.05 ±0.7 


0.405 ±0.16 


[27] 


777T + 7T _ 


0.148 ± 0.003 ±0.03 


0.085 ±0.07 


[27] 


r)ir + tt~ ir + ir~ 


1.11 ±0.01 ±0.3 


0.58 ±0.16 


[27] 


K+R-tt" 


0.0283 ±0.001 ±0.0005 


0.08 ±0.031 


[29] 


K+K-ir°TT° 


0.142 ±0.0006 ±0.01 


0.25 ±0.056 


[28] 


K+K-ir+n- 


0.708 ±0.002 ±0.1 


0.845 ± 0.09 


[28] 


K+K-ir+n-TT 


1.32 ±0.005 ±0.2 


0.915 ±0.16 


[27] 


T)4> 


0.013 ±0.0002 ±0 


0.12 ±0.06 


[29] 




0.129 ±0.001 ±0.04 


0.105 ±0.066 


[27] 


K + K-K + K~ 


0.0442 ± 0.0001 ± 0.01 


0.086 ±0.018 


[28] 


AA 


0.0501 ± 0.0002 ±0 


0.0355 ±0.0125 


[30] 


AS 


0.0231 ± 0.0001 ±0 


0.0071 ± 0.0082 


[30] 


EE 


0.0186 ± 6e - 05 ± 


0.023 ±0.011 


[30] 


Not included in the fit 


' + - 

f] TT 7T 


0.073 ± 0.001 ± 0.02 


0.035 ± 0.0325* 


[27] 




0.199 ±0.001 ±0.03 


0.065 ± 0.02* 


[28] 


/r *(892)/S-7r + 


0.141 ±0.002 ±0.005 


0.29 ±0.032* 


[28] 


/l7r + 7T" 


0.24 ±0.001 ±0.035 


0.3 ±0.12* 


[27] 




0.263 ± 0.0009 ± 0.09 


0.2 ±0.05* 


[27] 


0/0 


0.0439 ± 0.0002 ± 


0.065 ± 0.03* 


[28] 



* Errors are statistical only. 



TABLE V. Comparison between fitted and experimental cross sections of exclusive multi-hadronic channels in e + e~ collisions 
at y/s = 2.6 GeV. The errors on theoretical cross sections oshm are the Monte-Carlo integration statistical error and the 
error owing to uncertainty on branching ratios of resonances contributing to the channel, respectively. The experimental cross 
section values have been obtained by averaging available measurements, see text. 



